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We show that a vibrational instability of the spectrum of weakly interacting quasi-local harmonic 
modes creates the maximum in the inelastic scattering intensity in glasses, the Boson peak. The 
instability, limited by anharmonicity, causes a complete reconstruction of the vibrational density of 
states (DOS) below some frequency cuc, proportional to the strength of interaction. The DOS of 
the new harmonic modes is independent of the actual value of the anharmonicity. It is a universal 
function of frequency depending on a single parameter - the Boson peak frequency, ui, which is 
a function of interaction strength. The excess of the DOS over the Debye value is oc oj* at low 
frequencies and linear in lj in the interval cji <C a; <C Wc. Our results are in an excellent agreement 
with recent experimental studies. 

PACS numbers: 61.43.Fs, 63.50+x, 78.30.Ly 



I. INTRODUCTION 

One of the most striking properties of glasses is a maxi- 
mum in the inelastic scattering intensity observed in neu- 
tron and Raman scattering at frequencies between 0.5 
and 2 THz, far below the Debye frequency [1]. This so 
called Boson peak (BP) indicates an excess of low fre- 
quency vibrations over the Debye value which is given 
by the sound waves. It is seen as a maximum in the 
vibrational DOS divided by uj^, G{uj)/llj'^ but not neces- 
sarily in G{uj) itself [2]. The BP is a general feature, its 
magnitude varies strongly between materials. 

Despite numerous efforts, the BP remains one of the 
most intriguing problems of solid state physics. Some au- 
thors attribute the BP to vibrations of clusters of atoms 
of typical sizes [3,4]. The physical origin of these clusters 
in homogeneous amorphous media remains unclear and 
they have not been identified in numerical simulations. 

Another popular qualitative explanation of the Boson 
peak is a softening of acoustic phonons by static disorder 
[5,6] due to elastic Rayleigh scattering. However, even 
the most optimistic estimates show that the Rayleigh 
contribution is at least 4 times too small to explain the 
experimental data on thermal conductivity in glasses [5] . 
This mechanism is also in contradiction to the linear dis- 
persion law for acoustical phonons at the Boson peak 
frequency seen in molecular dynamics in Ref. [7] . 

Sometimes the BP is related to low lying optic modes 
of parental crystals [8-10] . Whereas the BP is a general 
feature of glasses such crystal structures with soft optic 
modes cannot be identified always. Such a mechanism is 
possible in some cases. Yet it remains a puzzle how these 
crystalline peaks are transformed in the glassy state to a 



shoulder in the vibrational DOS. 

Recent work on harmonic lattice models demonstrated 
that softening of disordered force constants can smear 
and push to low frequencies peaks which exist in the 
crystalline DOS [11-13]. In another approach the vibra- 
tions of a random distribution of atoms, interacting with 
a Gaussian-shaped pair potential, was studied [14] in a 
harmonic scalar approximation. Reducing the density 
the system becomes unstable. Approaching this insta- 
bility a low frequency peak appears in G{uj)/uj'^ which 
resembles the BP. The main drawback of these models is 
their neglect of the static displacement of the atoms in 
response to disorder. In real glasses where short range 
order is conserved such feedback always occurs since the 
forces between the atoms strongly depend on their dis- 
tance. In particular, the above models have no built in 
mechanism to stabilize vibrations with negative w^, un- 
stable modes. 

The proposed models of the Boson peak do not ac- 
count for anharmonicity effects which, as we will show 
in the present paper, become very important especially 
for small force constants. Glasses, at low frequen- 
cies/temperatures, are highly anharmonic as seen in most 
of their macroscopic thermodynamic functions. Anhar- 
monicity and static displacements, together, stabilize 
otherwise unstable vibrational modes. This does, how- 
ever, not imply that the vibrations at the Boson peak are 
anharmonic. Anharmonicity is essential in forming the 
equilibrium structure, and thus the force constants which 
determine the proper harmonic spectrum of low frequency 
modes. 

Another important point is that the previous explana- 
tions do not relate the Boson peak, which is one of the 
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universal properties of glasses, to other universal prop- 
erties, such as the two-level systems which dominate the 
low-temperature behavior or the plateau in the thermal 
conductivity at moderate temperatures. This relation- 
ship naturally emerges in our approach. 

II. QUASI-LOCALIZED VIBRATIONS (QLV) 

We present a new universal mechanism for the forma- 
tion of a BP in glasses, out of an originally flat DOS. 
This arises from the following 3 features: 1) quasi-local 
vibrations (resonant states) (QLV) with a smooth, struc- 
tureless initial DOS, (?o(w) at low lu, 2) an elastic interac- 
tion between them and 3) stabilization by anharmonicity 
when the system becomes unstable because of interac- 
tion. 

As the well known two level systems, the QLV are a 
typical feature of disordered systems. They arc addi- 
tional modes and are characterized by a large vibrational 
amplitude of some group of atoms. Their existence in 
glasses was predicted in Refs. [15-17]. They can be de- 
scribed as low frequency harmonic oscillators (HO) which 
couple bilinearly to the sound waves, see Refs. [18-20]. 
This in turn leads to a dipole-dipole interaction between 
different HO. The importance of the elastic interaction 
between local defects in glasses has been stressed by Yu 
and Leggett [21] and Grannan et al. [22]. 

The microscopic origin of QLV in disordered systems 
varies and depends on the type of disorder. Broadly 
speaking, the QLV can be divided into two groups. First 
there are materials where QLV exist independently of 
the structural disorder typical for structural glasses and 
amorphous solids. Orientationally disordered ("plastic") 
crystals belong to this group [23-25]. In these materi- 
als some molecular groups librate with low frequencies. 
In harmonic approximation these soft librations can be 
identified with QLV. The local potentials for the libra- 
tional motion can vary from site to site. In this case we 
have a distribution of the librational frequencies. The li- 
brations couple to the sound waves which in turn induces 
an interaction between them [22]. It depends on details 
of the material whether the interaction is strong enough 
to reconstruct the original spectrum of QLV completely 
or only partially. This effect is seen in recent measure- 
ments of dielectrtric loss spectra of ortho-carborane [26] . 
Similarly coordination defects in covalent materials can 
lead to QLV as was observed by Biswas et al. [27] in 
a simulation of amorphous Si (see also recent works [28] 
and [29]). QLV can originate from numerous defects such 
as off center ions or interstitial atoms. Depending on the 
"size" of the defect the QLV involve more or fewer atoms. 
Interstitial atoms are the prototype of a topological point 
defect. QLV of interstitial atoms in FCC metals were 
studied extensively in the past, see Ref. [30] for a review. 
These QLV have effective masses of four atomic masses 



and the crystal structure is strongly distorted by the de- 
fect and the low frequency of the librational QLV can be 
traced to the local strain. Already low concentrations of 
these interstitials are sufficient to destroy the crystalline 
structure completely. This is utilized in the interstitialcy 
model of glass formation [31]. 

This leads directly to the second group where the QLV 
result directly from disorder. Such modes are regularly 
found in computer simulations, e.g. for soft spheres [32], 
SiOz [33], Se [34], Ni-Zr [35], Pd-Si and Au-Si [36], 
NiB [37], in amorphous ice [38], in amorphous and qua- 
sicrystalline Al-Zn-Mg [39] and in simple dense fluids [40]. 

In these simulations the QLV were observed as local- 
ized vibrations with frequencies below the minimal sound 
wave frequency allowed by the size of the simulated sam- 
ple. A simple indicator of these "size-localized" QLV is 
the scaling of the participation ratio with system size. 
Increasing the system size the minimal sound wave fre- 
quency drops and the QLV are no longer localized in the 
simulation but show the typical properties of resonant 
modes, i.e. of low frequency local vibrations which cou- 
ple bilinearly to the sound waves. The exact eigenvectors 
of the interacting system of QLV and sound waves are su- 
perpositions of these two types of modes. The effect of 
system size on the appearance of QLV in simulations was 
discussed in detail in Ref. [41]. There it was shown that 
the exact eigenvectors at frequencies up to and above 
the boson peak can be decomposed into extended sound 
wave like modes and the local cores of QLV. The latter 
correspond to the harmonic oscillators of this work. 

The physical origin of these disorder induced QLV can 
be traced to local irregularities of the amorphous struc- 
ture. In dense packed metallic glasses these originate 
e.g. from the conflict of the local dense packing (icosahe- 
dral packing) and global dense packing (fee or hep) [42]. 
These local irregularities can be seen as centers of local 
strains [43]. 

Wc expect such local strain centers to be ubiquitous. 
The strains will have broad distributions, which will lead 
to broad distributions of QLV frequencies. Whereas lo- 
cal strains and QLV will be a general property of glasses 
the atomistic structure of QLV reflects the structure of 
the considered material. In dense packed metallic sys- 
tems the cores of the QLV have been found to be chain- 
like [44,41]. In Si02 they arc formed by a coupled rota- 
tions of Si04 tetrahedra [17,45]. In Se one has coupled 
chains and rings [34], etc. 

Another possible (and natural) mechanism of the QLV 
formation are low lying optical modes in parental crys- 
tals. Disorder in amorphous material would destroy the 
long range coherence of optical modes. This makes them 
practically indistinguishable from quasi-local modes. 

Together with the tunneling systems, the QLV form 
the main ingredient of the sofl potential model [46,47] 
(see review [48]). They manifest themselves in experi- 
mental values, e. g. the excess specific heat [47,49-51] 
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and the plateau in the thermal conductivity [20,52], in 
inelastic light [53] and neutron scattering [17] and they 
are observed in many numerical calculations, cited above. 



III. DENSITY OF STATES OF 
QUASI-LOCALIZED VIBRATIONS (QLV) 

One might think a priori that the QLV can have an 
arbitrary DOS, g{uj), depending on the particulars of the 
glass. We will show, that due to the interaction between 
the HO, g{w) is a universal function at low frequencies. 
This universality stems from the vibrational instability of 
the spectrum which occurs in nearly all systems of inter- 
acting HO. Anharmonicity stabilizes the system in new 
minima, and thus reconstructs the DOS to a new har- 
monic spectrum. As a result g(oj)/ijj'^ acquires a maxi- 
mum without a peak in g{w) itself. 




FIG. 1. The Boson peak, Eq. (1). 

Below we will derive the following form of the reduced 
DOS of these harmonic resonant modes (excluding the 
Debye part at low to which is not seen in Raman scatter- 
ing) 
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where C is a constant (Eq. (22)) and 



Z2{l^) 



^i,2(w) = -V V9 + 8(wVw)6±3. 
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(2) 



The function g{Ltj)/iJ^ is plotted in Fig. 1. It depends 
on a single parameter, w* characterizing the position of 
the Boson peak. The maximum of g{u))/u)'^, the Boson 
peak, is at loi, « l.lw*. For small frequencies, lo ^ Ub, 
g{uj) oc uj'^ while for large ones, uj ^ Wb, g{uj) cx uj. 



IV. VIBRATIONAL INSTABILITY 

To illustrate our central idea of a vibrational instabil- 
ity, wc start with a pair of interacting HO immersed in 
an elastic continuum. The potential energy is given by 

C4ar(a;i,a;2) = Mia;?a;?/2 + M-^Jixlj^ - hixxxi. (3) 

Here xi_2 are the HO coordinates, M\^2 the masses and 
a;i,2 the bare frequencies of the two HO, i. e. neglecting 
the bilinear interaction. The interaction strength is given 
by [53] 



'12 



(4) 



where gi2 accounts for the relative orientation of the HO, 
ri2 is their distance, p is the mass density of the glass 
and V is a sound velocity. The interaction between the 
HO is due to the coupling between a single HO and the 
surrounding clastic medium (the glass). This HO-phonon 
coupling has the form [20] 

Wint = Aa;£, 

where A is the coupling constant and e the strain. 
Diagonalization of Eq. (3) yields two frequencies 



Wl,2 




M1M2 



(5) 



The smaller value, , becomes negative when the inter- 
action /12 = / exceeds the threshold (critical) value 



7c = wia;2-\/MiM2. 



(6) 



This instability persists also in a system of many inter- 
acting HO. In a real physical system, anharmonic forces 
always stabilize an embedded HO in a nearby minimum 
of the potential energy. The position of this minimum 
depends on the interaction between HO. We are thus 
confronted with the many-body problem of finding the 
minima of the potential energy for a system of interact- 
ing anharmonic oscillators, similar to the one considered 
in Ref. [22,54]. The new frequencies, in these new min- 
ima, are real and different from the original ones. The 
harmonic vibrational spectrum is reconstructed. We will 
call this anharmonicity limited vibrational instability. 



V. STABILIZATION BY ANHARMONICITY 

Wc will now show that for weak interaction, /, the re- 
constructed DOS has, below a characteristic frequency 

oc |J|,a universal form irrespective of its original form. 
First, due to interaction, it becomes a linear function of 
frequency, g{w) oc w. Secondly, the displacements of the 
previously unstable oscillators from their old equilibrium 
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positions create static random forces which cause a sec- 
ond reconstruction of the DOS below another frequency 
cub ^ i^c- Due to so called sea-gull singularity [47] at 
Lu = the linear DOS is reconstructed to g(u}) cx oj'^ for 
u) << LVb- Together, these two reconstructions produce a 
maximum of g{u})/ui^ at a; = cOb- 

Let us consider a number of randomly distributed, in- 
teracting HO with concentration no and an initial DOS, 
5o('*^) (normalized to unity), in the frequency range from 
to ojQ, where goiuj) is a monotonously increasing func- 
tion of oj. For the harmonic part of the interaction we 
take the generalization of Eq. (3) and add an anharmonic 
term to stabilize the system 



C/anhar = (1/4) J2 ^i^t ^ > 0- 



(7) 



We will take the interaction I to be the small parame- 
ter of our theory, i.e. we assume that the typical random 
interaction, /, between neighboring HO is much smaller 
than the typical values of MlJq. As \I\ -C MlOq fre- 
quencies of order wq will be practically unaffected by the 
interaction whereas HO with frequencies w < will be 
displaced to new minima, where 



Wc ^ \I\/Mujo «; Wo- 



(8) 



Since the concentration of unstable HO is much smaller 
than the one of the stable ones a low frequency oscillator 
is typically surrounded by high frequency ones. We can 
simplify our consideration by again considering pairs of 
HO, one with a low frequency Wi < Wc and the other 
one from the bulk of the HO with a frequency a;2 of the 
order of cjq (sec Appendix A for the general case). Due 
to the combined action of interaction and anharmonicity 
the two HO will be displaced into new minima, xxq and 
X2O1 given by the equations (/12 = /) 



1x20 = xw{Miuf + Aixlo), 

Ixio = X2o{M2iol-\- A2xIq). 



(9) 



For |/| > Ic we need the nonzero solutions of these 
equations. Expanding around either minimum we find 
the new (harmonic) frequencies from the secular equation 



ai - Miw^ -I 

-I Q!2 - M20j'^ 







with 



ai = Miujf + 3A,x^Q, i = l,2. 



(10) 



(11) 



^Prom the condition ivi <C W2 follows 0:20 ^ Xio and, 
therefore, the term A2X20 in Eq. (9) can be neglected 
giving 



X20 = {I/M2UjI)xio. 



(12) 



a;io = a;iv/Mi7I7V(W-l- 



(13) 



As a result wc get from (10) under the condition 
'jJi/lo2 ^ 1 with this accuracy the new frequencies 
u}2 = i^2 and 



= 2a;2 [(//4)2 - 1] . 



(14) 



The smaller frequency (14) is the solution of the linear 
equation (compare with Eq. (10)) 



M2a;|(ai - Miw^) = 7^. 



(15) 



It is remarkable that for weak interaction the strength 
of the anharmonicity Ai does not enter the renormalized 
frequency (14). 

Near the threshold where (|/| —Ic)/Ic ^ 1, the smaller 
frequency squared is proportional to (|/|— Jc)//c. Pro- 
vided the distribution of the random quantity I is smooth 
one gets, therefore, below ujc a linear DOS {g{i-^) oc lu) ir- 
respective to the initial form of 5'o('^)- In Appendix A it 
is shown that the same result holds if one has a low fre- 
quency HO surrounded by several high frequency ones. 
Our numerical calculations (see Section VH) also show 
that this case is typical. 



VI. THE BOSON PEAK 

If the low-frequency HO with their reconstructed linear 
DOS were isolated, the problem would be solved. There 
is, however, a further interaction between these oscilla- 
tors which we have not taken into account so far. The 
low- frequency HO, displaced from their equilibrium posi- 
tions, create random static forces /. The force /, exerted 
on the ith oscillator by the jth one is 



fi — lijXjQ, 



(16) 



In a purely harmonic case, these linear forces would not 
affect the frequencies. Anharmonicity, however, renor- 
malizes the low frequency part of the spectrum [47], a 
manifestation of the sea-gull singularity. Consider an an- 
harmonic oscillator under the action of a random static 
force / 



U{x) = Ax^/A + MujIx^/2 - fx 



(17) 



where loi is the oscillator frequency in the harmonic ap- 
proximation. The force / shifts the equilibrium position 
from a; = to xo 7^ 0, given by 

Axl + Mu^xo - / = 0, (18) 

where the oscillator has a new (harmonic) frequency 

ojI,^=wI+3AxI/M. (19) 



and 
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If (wi) is the distribution function of frequencies loi 
and P{f) is the distribution of random forces, then the 
renormaUzed DOS is given by 



5M = J 9i{^i)d^i J dfP{f)6{: 



U) — LJ„ 



(20) 



As the forces between the HO are proportional to 



their sum is Lorentzian distributed (see Appendix B): 

(21) 



Assuming a; <C Wc and integrating Eq. (20) with gi{toi) = 
Ctoi we arrive at the integral 



g{Lu) _ 6C / cj \2 



dt 



with 



(22) 



(23) 



and after integration finally Eq. (1) is obtained. 

For small frequencies, below the Boson peak, lo <^ 
only small random forces / contribute to the second in- 
tegral in Eq. (20). In this case the distribution function 
P{f) can be approximated by a constant value, P(0), 
and we get from (20) 



(?(a;) (X lu^ dull 



UJi 



(24) 



As a result at low frequencies the renormalized DOS, 
g{ijj) cx uj^, Ref. [47]. For sufficiently large frequencies, 
LO LUb the action of random static forces on the HO 
spectrum can be discarded. In this case we recover the 
linear DOS, g{uj) cx uj. 

The BP frequency, LOb ~ l-lw*, is determined by 
the characteristic value of the random static force Sf, 
Eq. (16), acting on an HO with the characteristic fre- 
quency u)c- According to Eq. (16), it is due to the in- 
teraction between HO with frequencies of order of Wc, 
i.e. 



Ac) 

ij 



Juc, J « //no « Mujc/nogo{uJo) 



where Uc ~ no5'o(wc)i^c is the concentration of these HO. 
The characteristic displacement of a low-frequency HO 
(13) from the equilibrium position is XjQ ~ ujc^/ M/A. 
As a result, we get the estimate 



(25) 




5'o('*^o) 



and according to (23) 



1 /3 

UJb « [5o(<^c)/fl'o(<^o)] , Wfe < UJc- 



(26) 



Again, in lowest order the anharmonicity A does not en- 
ter this formula. 

As a result we get a following estimate for the recon- 
structed DOS 

{5o(w), uj > UJc, 

uj go{u!c)/uic, < w < Wc, (27) 
^^(Wo)/^^, UJ <UJb- 

If the DOS of the noninteracting oscillators is given 
by a power law, 50 (w) oc w", the BP frequency Ub scales 
with the interaction strength / as 



UJb 



(28) 



Since in accordance with (8) u)c oc |/|, we have always 
u)b <C UJc for n > 0. 



VII. NUMERICAL SIMULATION 




FIG. 2. Simulated density of states (.go(t^) oc uj'^, 
N = 2097) in a log-log representation. The arrows indicate 
the two characteristic frequencies uJb and uJc- 

To test our ideas by numerical simulations, we placed 
A'' oscillators with frequencies < < 1 on a simple 
cubic lattice with lattice constant a — I and periodic 
boundary conditions. To simulate random orientations 
of the oscillators we took for g^j, Eq. (4), random num- 
bers in the interval [—0.5,0.5]. The masses. Mi, and an- 
harmonicity parameters, Ai, were put to 1. The DOS for 
the noninteracting oscillators was taken as go{u)) oc w", 
with n = 1, 2, 3. 

Using the potential energy given by the generalization 
of Eq. (3) plus the anharmonicity (Eq. (7)) we then min- 
imized the potential energy, and in the usual harmonic 
expansion around this minimum calculated the DOS for 
different interaction strengths, J. This was repeated for 
up to 10000 representations. To check for size depen- 
dence we did the calculations for different N. Apart from 
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the case J = 0.07 the results did not change between 
N = 2097 and N = 4096. 

The predicted change-over in the w-dependence of the 
DOS at two characteristic frequencies lUc and loi, and the 
hnear part in between can be clearly observed in a log- 
log representation, Fig. 2, for 5(o(w) oc u)^ and J = 0.1. 
^From the calculated eigenvectors we find that, as ex- 
pected, at the lowest frequencies the HO are weakly cou- 
pled whereas near and above ui, the eigenmodes are com- 
plicated superpositions of many HO. 




FIG. 3. Simulated g{uj)/Lo^ for different interaction 
strengths {go{uj) ocoj^, N = 2097 and N = 4096 (J = 0.07)). 
The insert shows the scaHng of the crossover frequencies u>c 
(o) and iOb (•) with interaction strength J. 

Fig. 3 shows the dependence of the simulated g{u})/u}^ 
on the interaction strength J. We can see the general 
increase of W6 and related decrease of the BP intensity 
with increasing J. Our simulations cover one decade in 
BP frequencies. The insert shows that, in full agree- 
ment with our predictions (see Eq. (8) and Eq. (28)), the 
crossover frequencies change with interaction J as Wc oc J 
and Ub oc J^+^Z^. 
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FIG. 4. Boson peak in reduced units: Eq. (1) (solid line) 
and Raman data for lithium borate glasses [55]. The positions 
of the Boson peak (for different compositions x) are given in 
brackets. 
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VIII. DISCUSSION AND COMPARISON WITH 
EXPERIMENT 

In Fig. 4 we compare our theoretical curve, Eq. (1), 
with Raman scattering data of lithium borate glasses [55] 
with different compositions. The agreement is remark- 
ably good over the whole composition range. This sup- 
ports the idea of a miiversal shape of the Boson peak [56] . 
The shift of the BP to higher frequencies with increasing 
concentration of Li20 can be explained by an increase of 
the total concentration of QLV and consequently of their 
interaction. 

The agreement between theory and experiment is not 

confined to this class of material. This is exemplified 
by Fig. 5 which shows a comparison of the theoretical 



curve with numerical simulation results and neutron and 
Raman scattering data for 3 different glasses. 

One of the most important results of our theory is the 
predicted linear frequency dependence of the density of 
vibrational states above the Boson peak. It stems from 
the vibrational instability of interacting harmonic modes. 
Such linear behavior has been observed in many numeri- 
cal simulations on different glasses and model disordered 
systems [32,35,36,39,58,59]. It is also in a good agreement 
with many experimental results [60,57,61-66] where the 
vibrational DOS has a section with near linear frequency 
dependence. Fig. 6 shows this for vitreous silica at dif- 
ferent temperatures. Above the Boson peak the DOS 
increases approximately linearly with frequency. 

In this paper we dealt with the case of weak interaction 
between HO. If the interaction is increased the charac- 



6 



10-2 


J 1 m 

^ e H, 

p o ° n^^o^, 

■ n ^ 

: 


1 J 

■ 51 K 

O 1 55 K j 
□ 1 673 K 1 






Ok 




i a-SI02 
' 


ok \ 

o ; 

^ 




10° 






u (THz) 





FIG. 6. Density of states of vitreous silica for three tem- 
peratures taken from inelastic neutron scattering data [57]. 
The slope of the straight line on the figure is equal to -1.1. 

teristic frequencies Wf, and lOc grow and the gap between 
them narrows and finally disappears. Then our "BP" 
in g{uj)/uj'^ superimposes the "boundary peak" in 50 ('-^) 
at the edge of the assumed spectrum of naked (nonin- 
teracting) QLV. The BP can no longer be distinguished 
from the boundary peak or from a possible equivalent 
maximum in gii{Lo). It is possible that in some cases, 
e.g. orientational glasses, go{Lo) has a pronounced peak 
which is still visible after the reshaping of the DOS by 
interaction. 

Similar models with strong coupling between oscilla- 
tors were previously investigated by molecular dynamic 
simulations [22] and by the replica method [54]. For ex- 
ample in the Rcf. [22] Grannan et al. assumed that the 
dynamics "is completely dominated by the interaction 
between defects mediated by the strain field" . The local 
field (harmonic and anharmonic contribution) was ne- 
glected. We consider the opposite case where the local 
field (harmonic and anharmonic) dominates and interac- 
tion between dipoles is weak. 

IX. CONCLUSION 

In conclusion we present a universal picture of the 
BP formation in glasses. We have shown that the low 
frequency quasi-localized harmonic modes in glasses are 
destabilized by the weak bilinear interaction between 
them. Anharmonicity stabilizes the system in a new min- 
imum of configuration space. It completely reconstructs 
the low frequency part of the spectrum (at uj < ujc) and 
the Boson peak feature (at Ub <^oJc) naturally emerges. 
The thus created boson peak has a material indepen- 
dent shape. At low frequencies, below the BP, the vibra- 
tional DOS increases as g{Lo) oc co'^, and above the BP, as 
g{(jj) oc (J. 

Although the anharmonicity is responsible for this ef- 



fect, the final spectrum of stable vibrations remains har- 
monic. A remarkable feature of the presented theory is 
that the strength of the anharmonicity does not enter the 
new stable spectrum at all. It looks as if the anharmonic- 
ity does all the work, it stabilizes the system in a new 
minimum and reconstructs the spectrum and then, like a 
ghost, disappears. Therefore, the discussed phenomenon 
is independent of the variation of the anharmonicity be- 
tween different materials. The only parameter entering 
the final density of states is the strength of interaction 
between the HO. 

Compared to previous work, a main new result of our 
approach is the natural emergence of the BP on the un- 
structured, flat low frequency part of the initial spec- 
trum go{u))/w'^ where the DOS previously had no peaks. 
For small interactions the BP frequency is much smaller 
than the Debye frequency value. It shifts with interac- 
tion strength, 7, which explains the large variety of BP 
magnitudes found in experiment. In contrast with pre- 
vious models, a natural connection of the Boson peak 
phenomena with other universal properties of glasses is 
established. 

The authors VLG and DAP gratefully acknowledge the 
financial support of the German Ministry of Science and 
Technology and the hospitality of the Forschungszentrum 
Jiilich where part of the work was done. 

APPENDIX A: CLUSTER APPROACH 

Consider a cluster containing a low frequency oscillator 
with frequency lui < lUc surrounded by a large number, 
s — 1, of HO with much higher frequencies u)j ^ loq. 
Inclusive their interaction, the total potential energy of 
the cluster is 

Utot = }2 + 4 2] ^^^i - 2 • 

The equilibrium positions of the HO, xm, are given by 
the system of s nonlinear equations 

MiUiXio + AiXiQ = ^IijXjo, i = l,2,...,s. (A2) 

In the case of instability {xio ^ 0), in analogy to the 
previously considered case of a pair of oscillators, the 
static displacements of the high-frequency oscillators are 
much smaller than the one of the low- frequency oscillator, 
xio- Therefore, in leading order 

Xio = {hi /Miujf )xiQ, ij^l. (A3) 

Inserting these values into Eq. (A2) for i = 1 we get 

Miwfxw + Aixlo = xioY^j^- (A4) 
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Under the condition 



Miuj\ < k where k = 



(A5) 



the cluster becomes unstable and the low frequency os- 
cillator is displaced to a now minimum 



(A6) 



In the opposite case, Miajf > k, the cluster is stable and 

XiQ = 0. 

The eigenfrequencies of the interacting oscillators are 
the solutions of the secular equation of order of s 



-I2I OL-i — M'iiaP' 



-Is2 



Here the are given by Eq 



-lu 

-he 

as - M.u? 



:0. (A7) 



(11). In leading order in 



lij /Mujq the secular equation is a linear equation for uj 

(compare Eq. (15)) 

(ai - Mia;2) [] M,a,| - \{ M^. = (A8) 



or 



(A9) 



and the new low frequency of the system of coupled os- 
cillators is given by 



■^{Miojl-k), k<Miujl 
2 

— (fc-MiWi), k>MiJl. 



(AlO) 



As in the case of a pair of oscillators, the anharmonicity 
has been used in the derivation of Eq. (AlO) but does 
not appear in this or our final result, Eq. (1). 

To derive of the new (reconstructed) DOS the distri- 
bution of k, p{k), has to be calculated. Inserting Eqs. (4) 
and (A5) into the definition of p{k) gives 



p{k) 




(All) 



Here the angular brackets denote averaging over the po- 
sitions of the s — 1 high frequency HO, their frequencies 
and orientations. For simplicity we take equal masses 
Mj = M and for gij a uniform distribution in the inter- 
val [-1/2, 1/2]. 

Using the Holtsmark method [67] (see Appendix B) 
one gets 



where 



p{k) 



B = 



1 



B 



exp 



B^ 
2k 



TT Jhq I 1 



M. 



(A12) 



(A13) 



Here ng is the total concentration of HO in the cluster 
and (l/w)g is the w"^ moment of the normalized initial 
DOS gni'jj). This formula is a more accurate definition 
of the characteristic frequency lOc introduced in Eq. (8). 
Note that the distribution p{k) (A12) belongs to an im- 
portant class of one-sided stable distributions, Ref. [68]. 

Due to the combined action of interaction and anhar- 
monicity the DOS is reconstructed to §'(0;) = 2wG(lj^) 
with 

G(a;^) = (<5(a;^-^I3?)),^^^ 

y" dkp{k) j dtol Go{ujI)6 (w^ _ (ai4) 


and Go(ti^i) = go{u!i)/2u)i. Using Eq. (AlO) and inte- 
grating Eq. (A14) we obtain 



00 

G{uj^) = ^Jdkp(^k+ ^Mu;2^ Go 



00 

J dkpik)Go (^io^ + ^) . (A15) 


For low frequencies, a; ^ Wc , G(a;^) = const and g{u)) oc 
oj, i.e. the reconstructed DOS is a linear function of oj. 
For high frequencies the first term in Eq. (A15) can be 
discarded and the original DOS is reproduced, G(a;^) = 
Go(w^) for CO u)c- 

APPENDIX B: HOLTSMARK METHOD 

1. Distribution of random forces 

Let Xi be a random value with zero mean, (x) = 0, and 
finite {\x\) and let (TV ^ 00) be Poisson-distributed 
random points in 3-d space with concentration ric- The 
distribution function, P{f), of the random values 

i i 

can then be calculated by the Holtsmark method [67] . 
m =Mf-i:^))^^J dre^f^Fir) (B2) 
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with 



with 



F(t) = ( exp 



(B3) 



where angular brackets denote averaging over Xi and rj. 
Since the values ajj/r? are independent of each other 



F{t) 



-irx I 



(B4) 

which in the limit N ^ oo becomes 



F{t) = exp 



(B5) 



Using (x) = and changing the integration variable to 

y = \T\\x\/r^ we get 

Aim, 



F{t) — exp 



■M(N) I ^(i-cosy) 





• (B6) 



The integral equals tt/2 and Eq. (B2) is the Fourier trans- 
form of a Lorentzian distribution of random forces 



where the width of the distribution is given by 



(N>- 



(B7) 



(B8) 



2. p{k) distribution 

The same method can be applied to calculate the dis- 
tribution of the random quantity k 



E 



(B9) 



K{t) = ( exp 



(B12) 



Following the steps of the previous subsection we can 
write 



K{t) = exp 



-no 



(B13) 



and, introducing the new variable y = {x'^\T\)/{r^LO^), 
27r / \x 



K{t) = exp 



3 \ w 



r| (a + i/Ssignr) 

(B14) 



where 



oo oo 

/d f d 

-l^{l-cosy), 13 = j ^siny (B15) 



It is straightforward to show that a = j3 = v27r. There- 
fore 



with 



-ftr(r) = exp — i3-\/|T| (1 -h isignr) 



TT Fk / 1 

^=3V2^°^°W 



(B16) 



(B17) 



After integration in (Bll) with K{t) from (B16) we fi- 
nally get 



p{k) = 



2^^3/2 ^^Pl"^ 



(B18) 



where the random and uniformly distributed, 

—xq/2 < Xi < .To/2, Ti are {N oo) Poisson-distributed 
random points in 3-d space (concentration no) and Wj are 
random frequencies of HO distributed in the interval 



< oji < ojQ 



(BIO) 



with a DOS, go (t^) 7 normalized to unity. 

Analogously to (B2) the distribution p{k) can be writ- 
ten as 

oo 
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